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A quantized vortex in the Bose-Einstein condensation 
(BEC) , which is Icnown to be unstable intrinsically, is demon- 
strated theoretically to be stabilized by the finite tempera- 
ture effect. The mean-field calculation of Popov approxima- 
tion within the Bogoliubov theory is employed, giving rise 
to a self-consistent solution for BEC confined by a harmonic 
potential. Physical origin of this vortex stabilization is in- 
vestigated. An equivalent effect is also proved to be induced 
by an additional pinning potential at the vortex center pro- 
duced by a focused laser beam even at the lowest tempera- 
ture. The self-consistent solutions give detailed properties of 
a stable vortex, such as the spatial profiles of the condensate 
and non-condensate, the particle current density around the 
core, the whole excitation spectrum, and their temperature 
dependences. 

PACS numbers: 03.75.Fi, 05.30. Jp, 67.40.Vs 



I. INTRODUCTION 

Many-body dilute Bose systems are the focus in con- 
densed matter physics since the experimental realization 
of Bose-Einstein condensation (BEC) in alkali atom gases 
at ultra-low temperatures . There is a pile of recent 
theoretical and experimental works devoted to elucida- 
tion of BEC [QH]. Experimental works are mainly con- 
ducted for dilute Bose systems, whose atom number is 
typically ~ O(IO^) for ^■^Na and ^''Rb. Gases are con- 
fined magnetically in a harmonic trap. The BEC tran- 
sition temperatures are in a range of 0(/xK). The ex- 
perimental forefront is rapidly extending; one can tune 
the particle-particle interaction constant by utilizing the 
so-called Feshbach resonance phenomena Q, or one can 
confine a system purely optically so that several atomic 
hyperfine substates are simultaneously condensed, giving 
rise to an opportunity to explore a BEC with internal de- 
grees of freedom . Several proposals of experimental 
method to realize a vortex state have made theoretically 
|p|, pi)|JTll] . This paper concerns the vortex states in some- 
what different way. 

As for the earlier theoretical developments, the micro- 
scopic theoretical work on BEC has started with Bogoli- 
ubov [l^ ] long time ago and followed by several important 
progresses, such as Gross Q, Pitaevskii Q, lordan- 
skii [Hi, Fetter and Popov (l7|. These mean-field 
theories are particularly suited for treating the present 
gaseous BEC systems trapped optically or magnetically 
in a restricted geometry. Therefore, the current the- 



oretical studies pq l mainly focus on examining these 
mean field theories to extract the spatial structures of 
the condensate and non-condensate and low-lying collec- 
tive modes. Agreement between these mean field theories 
and experiments is fairly good 

Previously we have demonstrated that Bogoliubov ap- 
proximation (BA) of the mean-field theory cannot sustain 
a stable vortex This somewhat unexpected result, 

which had been hinted by DoddjES, is now clearly ex- 



plained physically by Rokhsar |21| who shows that the 
vortex is always unstable against the translational mode 
with zero energy. Thus the vortex core in a cylindrically- 
symmetric harmonic trap spirals out towards the outside 
of a system where the condensate density is low. 

The main purpose of this paper is to examine the sta- 
bility of a quantized vortex which is directly related to 
the question of superfluidity of the present BEC systems. 
This instability can be expressed numerically either by 
comparing the free energy of systems with and without a 
vortex (stability), or by checking whether an eigenstate 
with a negative eigenvalue exists or not in a system with 
vortices (metastability). The latter is far serious than 
the former and is discussed in this paper as "stability" . 

In the present paper, which follows our series of the 
papers |l9||2^ ] based on the same formalism and numeri- 
cal procedure, we show that the intrinsic instability of a 
quantized vortex can be remedied either by raising tem- 
perature (T) or by introducing the pinning potential at 
the center of a cylindrically symmetric rectilinear vor- 
tex line, which prevents the vortex from spiraling out. 
The latter method for the vortex stabilization is pointed 
out by several authors , whose experimental feasibility 
seems to be quite high since it may be easy to realize a 
desired pinning potential by sending a focused laser beam 
at the vortex center. Here we provide microscopic self- 
consistent calculations within the Popov approximation, 
yielding detailed properties of a stable vortex, beyond 
Gross-Pitaevskii theory which treats only the condensate 
component. This enables us to discuss at finite temper- 
atures. 

Various spatial structures of the quantized vortex in 
rotating BEC systems are investigated based on this 
scheme of stabilization. This investigation is important 
since the microscopic understanding of the vortex based 
on the self-consistent picture has been lacking so far, ex- 
cept for a few attempts p^-|2^ for present BEC systems. 
Several points are worth exploiting: (1) the spatial pro- 
files of the condensate and non-condensate, (2) the dis- 
persion relation of the excitation spectrum of the system. 
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and (3) the low-lying excitations localized around a vor- 
tex core, which is known as the Kelvin mode p6t . These 
features may be directly observable and should be com- 
pared with those in vortices in a superconductor where 
extensive studies have been compiled . 

Arrangement of this paper is as follows: In the next 
section we recapitulate the Popov approximation (PA) 
which is one of various types of mean-field theories (see 
for example, Refs. jl^j2^) due originally to Bogoliubov 
I p^ . We explain our numerical procedure to solve a self- 
consistent equation within the PA. In Sec. 3 the physical 
origin of stabilization is discussed and the detailed cal- 
culations for a stable vortex at finite temperatures are 
shown. We also compare the calculated vortex systems 
to the vortex state in type-II superconductors. The pin- 
ning potential is introduced to yield a stable vortex in 
Sec. 4 where the strength and width of the potential 
are critically examined for the vortex stabilization. We 
devote the final section to discussions and conclusion. 



II. FORMULATION AND NUMERICAL 
PROCEDURE 

A. Mean Field Approximation 

We start with the following Hamiltonian in which Bose 
particles interact with a two-body potential: 

H= dr¥{r) i — ^ + ^(^)"'^) 

-f| Jdr¥{r)i'\r)i/{r)i'{r) (1) 

where the chemical potential /i is introduced to fix the 
particle number and V(t) is the confining potential. The 
two-body interaction between an atom at ri and one at 
r2 is assumed to be gS{ri — with g being a positive 
(repulsive) constant proportional to the s-wave scattering 
length o, namely g = 47r?i^o/m (m is the particle mass). 

In order to describe the Bose condensation, we assume 
that the field operator ^ is decomposed into 



*(r) = 0(r) + V(r) 
where the ground state average is given by 
(#(r)) = 0(r). 



(2) 



(3) 



A c-number (f)(r) corresponds to the condensate wave 
function and ipir) is a q- number describing the non- 
condensate. Substituting the above decomposition (|^) 
in we obtain 



H : 



dr 



r(r)^/i(r)0(r) + -g|0(r)p 



(r) 



+V't(r) {h{r) + g|0(r)|2} 0(r) -t- h.c. 



-fV;+(r){/i(r) + 25|0(r)|2}^(r) 

+ |v3t(r)^t(r)^(r)^(r) + |^(r)^(r)^*(r)^*(r) 

+gij;\r)ip{r)ip{r)4)* (r) + gti;"^ {r)ip^ {r)ip{r)4){r) 
+ pt(r)^t(r)^(r)^(r)] (4) 



where h{r) 



V{r) — /i is the one-body Hamil- 



tonian. Let us introduce the variational parameter: the 
non-condensate density p(r) = (ip'^tp) ^^'^ approximate 
as tp^ip'^P ~ 2-0^, and ^'^^'^^ip ^ Aij^'^tpp. 
Then, H is rewritten as 



H = / rfr 



r(r)|Mr)<^(r) + ^ff|</'(r)p|</)(r) 

-KV't(r) {h{v) + 5|0(r)|2 + 2gp{v)) cj,{v) + h.c. 
+^\v) {h{v) -f 2.9(|0(r)|2 + p{v))] V-(r) 

+ pt(r)^t(i.)0(i.)0(r) +h.c. 



(5) 



In order to diagonalizc this Hamiltonian, the following 
Bogoliubov transformation is employed, namely, ■(/'(r) is 
written in terms of the creation and annihilation oper- 
ators rjq and 77J and the non-condensate wave functions 
Uq{Y) and Ug(r) as 



(6) 



where q denotes a set of the quantum numbers. This 
leads to the diagonalized form: H = Eq -I- J2q^qVl''h- 
The condition that the first order term in ^/'(r) vanish 
yields 



[h{r)+g\^{r)\^ + 2gpir)] 0(r) = 0. 



(7) 



When p(r) is made zero, it reduces to the commonly used 
Gross-Pitaevskii equation 



[h{r 



Mr)\'] 0(r)=O, 



(8) 



which is a non-linear Schrodinger type equation. 

The condition that the Hamiltonian be diagonalized 
gives rise to the following set of eigenvalue equations for 
Uq{r) and Vq{r) with the eigenvalue e^: 



[/i(r) + 2g{|0(r)|2+p(r)}] Uq{r) 

-g02(r)w,(r) = eqUq{r) 
[hir) + 2gmr)f+pir)}] Vqir) 

-g(t)*^{r)uq{r) = -eqVq{r). 



(9) 



(10) 



The chemical potential fi is contained in h of the left 
hand side of Eqs. (||) and ( p^ . Therefore p corresponds 
to the origin of the energy for each eigenvalue Eq. 

The eigenfunctions Uq{r) and Vq{r) must satisfy the 
normalization condition: 
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J ['u*p{r)uq{r) - v*{r)vq{r)] dr = Sp^q. (11) 

This also means that / juq^iir is larger than / |wqp(ir 
for an arbitrary q. The variational parameter p(r) is 
determined self-consistently by 

p(r) = (V^V^) 

= E [I^.WIV(£,) + Mr)\'{f{eq) + 1)] , (12) 

where /(e) is the Bose distribution function. Equations 
(0); (H); (0)' (0); ^^'^ dll) Constitute a complete set of 
the self-consistent equations for the Popov approximation 
(PA). The iterative calculations of these equations yield 
a convergent self-consistent solution. 

If the non-condensate component p{r) is neglected in 
this set of equations, the resulting set of the equations 
yields the so-called Bogoliubov approximation (BA). 
This was extensively discussed analytically by Pitaevskii 
0, lordanskii ||l||, and Fetter |l|]. 

The expectation value of the particle number density 
is given as 



(n(r)> = |0(r)|2+p(r), 



(13) 



that is, the total density consists of the condensate part 
(j){r) and the non-condensate part p(r). The particle cur- 
rent density is calculated as 

j(r) = :r^{0*(r)V^(r) - ^(r)V0*(r)} 

+ -^(^t(i.) . vV'(r) - VV;t(r) • ^(r)). (14) 

The local density of states is calculated as 

N{E, r) = E {l^^qMP + Sie^ ^ E). (15) 

q 

This quantity was observed directly by scanning tunnel- 
ing microscopy for a vortex in superconductors |p9| . 

B. Vortex Description 

We now consider a cylindrically symmetric system 
which is characterized by the radius R and the height 
L. We use the cylindrical coordinate: r = {r,9,z). We 
impose the boundary conditions that all the wave func- 
tions vanish at the wall r — R, which is taken far enough 
from the vortex center at r = 0, and the periodic bound- 
ary condition along the z-axis. When a vortex line passes 
through the center of the cylinder, the condensate wave 
function (f)(r) is expressed as 



where (j){r) is a real function and w is the winding num- 
ber, w = corresponds to the non- vortex case and w = 1 
to the vortex case. The w > 2 case is not considered here 
because this state is energetically too unstable. The non- 
condensate density p(r) is a real function, depending only 
on r, that is, p(r, 9, z) — p{r). It is also seen from Eqs. (||) 
and ( p^ that the phases of Uq(r) and Vq{r) are written 
as 



Wq(r) = Uq(r)e*«^^e*(«''+"')'' 



(17) 
(18) 



The set of the quantum numbers q in is described by 
{qr,qe, Qz) where qr = 1, 2, 3, • • •, qe = 0, ±1, ±2, •••,52 = 
0, ±2tt/L, ±Att/L, ■ ■ ■. Note that only u{r){qg = -1) and 
v{r)(qg = 1) are non- vanishing at r = 0. 

By following the method by Gygi and Schliiter [ pO| , 
the functions Uq{r) and Vq{r) are expanded in terms of 



Vyir) = . J\,A (aLi — 



Ji 



as 



i 



(19) 

(20) 
(21) 



where Ji, (r) is the Bessel function of i^-th order and a)> 
denotes «-th zero of J^. 

The eigenvalue problem of Eqs 
duced to diagonalizing the matrix: 



Mij(qe + w, g^) -B^j{q0,w) 
Bjjiqe, w) -A,j{qe - w, g^)y 



) and 


(0) 


are re 














\ ■■ ) 




\ ■ ) 



(22) 



for each qg and qz, where 



(16) 



+ / {V + 2gic^^+p))^i'^^y^rdr (23) 
Jo 

B,^,{qg,w)^ f gcf'^+^^^^^^dT. (24) 
The following symmetry relation should be noted: 



3 



(c, d, e) for qg <^4> (d, c, — e) for — qg. 



(25) 



We determined the signs of qg and e by the normalization 
condition Eq. (p^. 

It is not self-evident a priori that Eq. ( p2| ) gives real 
eigenvalues because the Hamiltonian matrix in Eq. ( |2^ ) is 
not symmetric. However, these eigenvalues can be proved 
to be real (positive or negative) by a way similar to that 
of Fetter |6|. 

The circulating current density jg (r) in Eq. (|lj) is ex- 
pressed as 



je{r) 






(26) 




h w 

m r 




(27) 




= Ay. 


{'^^ + "MOlV(s,) 






q 

<?e - 
r 


-KMP(/(eq) + i)} 


(28) 



where the total current consists of the condensate com- 
ponent i^f^\r) and the non-condensate j^\r). 



C. Calculated System 

We have performed self-consistent calculations of a gas 
of Na atoms trapped radially by a harmonic potential 
V{r) = im(27ri/r)^r^. We use the following parameters: 
m — 3.81 X 10~^^kg , the scattering strength a = 2.75nm, 
and the radial trapping frequency Vr — 200Hz. The 
area density per unit length along the z-axis is chosen 
to he Uz — 2 X 10*/ /im , so that the peak density will be 
O(10^°m~'^). The system size is set to R = 20^m and 
L = lO/Ltm. 

We have employed the energy cutoff method to calcu- 
late Eqs. (|l^) and (^8|) where the calculation terminates 
when the obtained eigenvalues exceeds about ^Ohv^. This 
method treats ^ 80000 eigenfunctions, although it is low- 
est two modes, namely the condensate state and the low- 
est core localized state (LCLS), that play significant role. 



III. STABLE QUANTIZED VORTEX AT FINITE 
TEMPERATURES 

The isolated single quantized vortex is unstable within 
certain mean-field approximations p9| , pl| in a non- 
rotating system. When the lowest excitation eigenvalue 
of the system becomes lower than the energy of the con- 
densate, the whole calculation based on the assumption 
that the energy for the condensate is lowest |l9j is inval- 
idated. This means that the vortex is unstable. 

In the finite T calculations shown below the quantized 
vortex becomes stable in the BEC systems confined by 
a harmonic potential above a certain temperature. We 



analyze how the vortex state at finite T can become sta- 
bilized, and examine fundamental properties of a stable 
vortex. 



A. Spatial Structures of Condensate and 
Non-condensate 

The spatial profiles of the condensate and non- 
condensate density for several temperatures are shown 
in Figs. |l|(a) and (b) respectively. It is seen from Fig. 
|l|(a) that the condensate |'^(?')p vanishes and changes 
quadratically (cx r^) near the vortex center r = 0, which 
is contrasted with the linear behavior (cx r) in the Gross- 
Pitaevskii theory Eq. (||). This unique shape of |(/)(?')J^ 
refiects the rounded shape of the 2gp{r) term in Eq. (|^). 
[This p{r) is depicted in Fig. ||(b).] 

This feature differs also from the superconducting case 
where the pair potential changes linearly |2^. As T in- 
creases, the core radius increases as seen from the inset. 
The condensate is pushed outwards and converted into 
the non-condensate. Thus the total number of the con- 
densate decreases monotonically as T increases. 

In Fig. |^(b) the spatial variation of the non-condensate 
density p{r) is displayed. Refiecting the core structure, 
the non-condensate fraction accumulates at the vortex 
center. The core radius grows as T increases as inspected 
from Fig. ^(b), which is analyzed in more detail shortly. 

The overall feature of the spatial change of the total 
density consisting of the condensate and non-condensate: 
|(/)(r)p -I- p{r) is shown in Fig. 0(c). As T increases, the 
spatial profile tends to uniformly distribute, and the max- 
imum position of the density from the core also increases, 
indicating that the effective core radius expands. The to- 
tal density does not vanish anywhere inside the confined 
potential even at the vortex core. It is also interesting 
to notice the existence of a focal point near the center 
which dose not move as T changes. 



B. Excitation Spectra 

In Figs, ^(a) and ^(b) we display the dispersion rela- 
tions of the eigenvalues along qg and qz respectively. The 
excitation modes are not symmetric about the g^-axis 
because of the presence of the circulation. However, at 
the higher energies the distribution of these eigenmodes 
become symmetric. The modes characterized hy qz ~ 
are shown as black dots in Fig. ^(a). The lowest eigen- 
mode aX qg — —I is very important in understanding the 
mechanism of the vortex stabilization. 

The modes with {qr,qg) = (1,-1) and various qz in- 
cluding the LCLS correspond to the so-called Kelvin 
mode known in classical vortex. In Fig. ^(b) we dis- 
play the dispersion relations of the eigenvalues with 
qg = —1 along g^. This dispersion relation is predicted 
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by Pitaevskii 1 14| in a spatially uniform system with den- 
sity riuni except near the core within the BA, at the long 
wave length limit as 



2„2 



e{qr = = -1,92 



In- 



1 



n q. 



(g^funi < 1), 



which is named the Kelvin mode, where 



n 



y/2mnunig 



(29) 



(30) 



The dispersion relation of the corresponding modes with 
qe = —1 is shown in Fig. ||(b) as the lowest edge. We 
have tried to fit this prediction Eq. (|29| ) with our result 
by adjusting nuni in Cuni of Eq. (|30|). The dotted curve 
in Fig. ^b) is when riuni = 0.1 x ripcak- The fitting is 
rather satisfactory. 

The spatial profiles of the eigcnfunctions with qg = —1 
are depicted in Fig. |[ There is no rapid spatial oscilla- 
tion in these eigenfunctions corresponding to the Friedel 
oscillation seen in the Fermi systems [^,^. It is also 
noticed that a slow oscillation of the condensate near the 
vortex core seen in several Monte Carlo computations |3^] 
for superfluid ^He is absent here. The present weakly in- 
teracting case differs from that in strongly interacting 
system in several respects. 



C. Vortex Stabilization at Finite T 

We discuss the origin of vortex stabilization. Let us 
start with the particle number density J(/)(r)p (conden- 
sate) and p(r) (non-condensate) in Fig. ||. It is the eigen- 
states that compose the whole system. The wavefunction 
(/)(r) with the energy /Lt corresponds to the condensate 
fraction and the non-condensate density p{r) is composed 
of the eigenstates with the eigenvalues and the wavefunc- 
tions Eq, Mq(r), and fq('r). The eigenstate with Sq = 
which equals corresponds to the condensate. 

The states with the angular momentum index qg — —1, 
whose relative motion to the circulating current of the 
condensate is at rest in the laboratory frame [qg + w in 
Eq. (|7|) is equal to zero], are likely to have the low- 
est eigenvalue among various states. It is the eigenstate 
q = {qr,qe,qz) — (1,-1,0) whose eigenvalue becomes 
negative and leads to the vortex instability. The spatial 
profiles of the eigenfunctions q = (1 to 5, —1, 0) have de- 
picted in Fig. ^. We can see from it that these states are 
localized near the core (r = 0). Black dots at qg = —1 in 
Fig. ^(a) correspond to the wavefunctions in Fig. ||. We 
can see that the lowest mode q = (—1,0, 1) is strongly 
localized at r = 0. Let us call this state the lowest core 
localized state (LCLS). The vortex state is (meta) stable 
as far as the LCLS exists with a positive eigenvalue. Now 



the question is why the energy of the LCLS is positive at 
finite T as shown in Fig. ||(a). 
We define the effective potential 



V,s{r) = V{r)+2g{\cj){r)f+p{r)), 



(31) 



whose combination appears in the eigenvalue equations 
Eqs. (|) and 

The eigenvalues of various modes are estimated by in- 
tegrating T4ff(r) with each of the wavefunction over the 
spatial volume. Of course the LCLS is not an exception. 
The shape of V^ff (r) near the core (r ^ 0) determines the 
eigenvalues of the LCLS. As seen from Fig. ^, V^ff(r) at 
the vortex center is pushed upward as T increases. This 
is why the LCLS has a positive eigenvalue. Although 
the vortex stability seems to increases with T, to deter- 
mine the stability conditions is a complex problem and is 
discussed in the next subsection. The increase of Vcsir) 
near the core is caused by the contribution from the non- 
condensate p{r). Thus Vcs{r), especially the term 2gp{r) 
m Eq. (|l]), works as an effective pinning potential which 
prevents the vortex from moving outwards. Then let us 
consider why the non-condensate p{r) has a peak at the 
core. 

The non-condensate density p(r) is composed of the 
wavefunctions u{r) and v{r) as seen from Eq. (|l|). Note 
that only u{r){qg = —1) and v{r){qg = 1) are non- 
vanishing at r = as mentioned before. The localiza- 
tion of p(r) at small r mainly comes from the LCLS 
term [q = (1,-1,0) term] in Eq. (|l|). Fig. | again 
shows that lit(?')l^ coefRcient of this term is strongly lo- 
calized at r = 0, and Fig. ^(a) shows that /(e) coef- 
ficient has huge value because the eigenvalue e is very 
small. Note that /(e) is the Bose distribution function 
and lim^^o /(s) = +oo. The LCLS, whose eigenvalue is 
approaching zero from above, itself forms a pinning po- 
tential to keep the eigenvalue positive combined with the 
Bose distribution function. We can see the wavefunc- 
tions and the eigenvalues in Fig. ^ which displays the 
local density of states defined by Eq. (|l5|). This may be 
helpful to understand the dominance of the LCLS in p(r) 
at small r. 

This mechanism of the T-dependent stabilization is 
based on the PA, which introduces the non-condensate 
density p(r) as an effective potential ignoring the friction 
between the condensate(w = 1) and the pinning part of 
the non-condensate [(90 -I- w) = 0]. In other words, the 
mechanism is based on superfluidity. 

It is also obvious that large enough Vcff(r) near the 
center is decisive in stabilizing the vortex, irrespective of 
the origin either coming from gp{r) or from V{r). The 
deformation of V{r) also stabilizes the vortex state and 
depend neither on the PA nor on the friction between 
p{r) and 4>{r). This deformation of V{r) is discussed in 
the next section, and the rest of this section is devoted to 
the various properties of the T-stabilized system which 
has been obtained. 
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D. Stability 



F. Local Density of States 



As mentioned before, to determine the stability con- 
ditions is a complex problem. It is natural to consider 
the eigenvalue £ > for LCLS and large enough Vcs{r) 
at r = as possible conditions for the stability. What 
makes the matter complicated is a relation of these two 
conditions. To express more explicitly, 



K=ff(0) =5|ULCLS(0)P/(£LCLS) 

+ other contributions from 
"u-ile = — l)'s and v{qe — l)'s. 



(32) 



The T-dependence of the selected lowest modes are dis- 
played in Fig. ^. As anticipated from the T-dependence 
of the lowest mode, which gradually decreases on lower- 
ing T, the whole calculation breaks down when this mode 
becomes negative from above. But this lowering process 
of e with T will not proceed smoothly because the Bose 
distribution function at small s could make the pinning 
part of Vcff(r) very large. 

When r = 0, the effective potential at the center 



(33) 



q('3e = l) 



is not likely to become large enough, although we cannot 
exclude its possibility completely. 

We only point out that these problems are open and 
shall not discuss them further in this paper. We failed 
to determine the self-consistent density profiles p{r) and 
|(/)(r)p for T < 200nK, which is the lowest temperature 
we can reach. 



E. Core Radius 

Figures |^(a), (b), and (c) show that the core radius 
increases as T increases. We define here the core radius 
by the distance ^ from the center where ^ max(|0(r)p) = 
I The T-dependence of this ^ is exhibited in Fig.||. 
The core radius ^{T) is seen to be T-dependent and de- 
crease as T is lowered. We define ^poak = ^/ y^2mnpcakff 
where ripoak = max(|(/)(r)p) is the peak density of the 
condensate. This resembles the classical definition 
of the coherence length ^uni in Eq. (|3^) for a uniform 
system. The core radius ^ is almost equal to 1.5 x ^poak 
as seen from Fig. ^ 

The total current density jg^^ (r) -I- jg^^ (r) defined in 



Eq. (26) is shown in inset of Fig. ||. The distance from 
the vortex center at which the current amplitude takes a 
maximum is another measure of the vortex core radius. 
This quantity shown in Fig. ^ also decreases upon 
lowering T. 



The local density of states N{r, E) given by Eq. ( p^ ) 
is shown in Fig. The states with qg ^ —1 are seen to 
localize at the vortex center, giving rise to several promi- 
nent peaks in the local density of states. 

These features around the core may be observable and 
should be compared with those in vortices in a super- 
conductor where extensive studies have been compiled 
pTf . This is particularly true for the localized excita- 
tions near a vortex core studied by Caroh et al. p5| , ^ . 
These core excitations are now recognized to be deci- 
sive in understanding the fundamental properties of the 
vortex. The Popov theory just corresponds to the so- 
called Bogoliubov-de Gennes theory in their mathemati- 
cal structures, which is widely used to treat various spa- 
tially non-uniform superconductors having interface or 
vortex lolJII]. 

We also note that in superconductors N{r, E) is di- 
rectly observed by scanning tunneling microscope |2g| ] 
and analyzed theoretically within the similar theoretical 
framework po| , ^j33| quite successfully. 

The characteristic core radius is estimated as the co- 
herence length ~ 0.5/im which is contrasted with a few 
A in superfluid ^He. Thus there is a good chance to in- 
vestigate the detailed core structure of the present BEG 
systems. 



G. Circulating Current Density 



,•(1) 



(2') 

) + jg (r) defined in 
The non-condensate 



The total current density jg 
Eq. ( p6| ) is shown in inset of Fig 

contribution jg (r) and its qg components are depicted 
in the main panel of Fig. H. In the immediate vicinity of 

(2) 

the core jg (r) is governed by the qg = and the qg — 
—2 components. The positive (negative) qg^s give rise 
to the positive (negative) contribution to jg{r), except 
for the qg = —\ which has positive contribution. The 
|u(r)p term is dominant in most qg of Eq. (p8|). This is 
why positive (negative) gg's give the positive (negative) 
contribution. Since we are treating systems with w — 1, 
the S£±2£ coefficients of the |M(r)p terms become zero 
and the |t'('')P terms with positive sign become relevant 
when qg = —1. This is why the qg = —1 component 
gives positive, not negative, contribution. The net result 
gives a finite contribution to superfluid flow from the non- 
condensate. 



IV. VORTEX STABILIZATION BY PINNING 
POTENTIAL AT T = 

It is now clear that introduction of a pinning potential 
at the vortex center stabilizes the quantized vortex even 
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at T=0. Here we explicitly demonstrate it by a model 
pinning potential described by 

l^pi„ning(r) = V^e'^o (34) 

in addition to the harmonic potential. This pinning po- 
tential may be realized by focusing a laser beam on the 
vortex center. Since we are imposing the periodic bound- 
ary condition along the z-axis, topologically the system 
is equivalent to a toroidal geometry if the system length 
L is large enough. 

We have done self-consistent calculations at T = 
under the presence of the pinning potential as the same 
method as in the previous section. The self-consistent so- 
lution for Vq = 5 trap unit and tq = l/im and the other 
parameters are same as before is depicted in Figs, ^(a), 
(b), and (c). It is seen from those that the vortex is 
stable because the LCLS is now all positive as shown 
in Fig. |^(b), allowing the fully self-consistent solution. 
We note that the eigenfunction u{r){qr = 1) shown in 
Fig. ||(c) has still large amplitude at r = 0. The conden- 
sate is pushed outwards because of the pinning potential, 
and its initial rise at the vortex core is steep [Fig. §(a)]. 
The non-condensate fraction is vanishingly small in this 
example at T = 0. 

In a larger pinning potential (Vq = 50 trap unit and 
tq = 1/im) exemplified by Fig. the condensate is far- 
ther pushed outwards and the vortex center is empty 
space [Fig. 0(a)]. It resembles a ring-shaped EEC sys- 
tems which has been considered theoretically |^ , ^ . 
The wave functions of the states q = (1 to 5,-1,0) in 
Fig. |To|(b) is squeezed at small r by the pinning poten- 
tial. There is no LCLS. This deformation of u{r) directly 
lifts the eigenvalue of the lowest eigenstate at — — 1 
in Fig. |o|(b) to stabilize further the vortex state. The 
overall distribution of the excitation spectrum is almost 
symmetric about the qg axis [Fig. |l^(b)]. All the corre- 
sponding eigenfunctions with qg = —1 vanish at r = 
[Fig. p^ (c)], which is contrasted with those in Fig. || and 
Fig. Ifc). 



realized by raising T or introducing an external pinning 
potential. In the former case the LCLS itself acts as the 
effective pinning potential, amplified through the Bose 
distribution function. The nature of the frictionless inter- 
action between the condensate and the non-condensate in 
the Popov theory enables this T-dependent mechanism. 

The detailed vortex properties are analyzed for the sta- 
ble vortex. Several eminent features are worth mention- 
ing: (1) The condensate |0(r)|^ behaves as oc near the 
vortex core and the non-condensate has a peak at the 
core. (2) The characteristic length or the core radius is 
well described by a formula h/ yJ2mn-pcakg with ripcak be- 
ing the peak density of the condensate. (3) The particle 
density is non- vanishing everywhere inside a system even 
at the core center when a system is T-stabilized. (4) The 
local density of states exhibits large peak structures at 
the core center for the low energy side, and the peaks 
mean the existence of the LCLS and other core localized 
states. These characteristics may be directly observed 
once the quantized vortex is realized in the present BEC 
systems. 

We hope that the vortex state in BEC of alkali atom 
gas should come true experimentally. If realized, we may 
have a chance to unify theories of the present dilute BEC 
systems and the strongly interacting superfluid ^He sys- 
tems. We also hope that the precise and detailed analysis 
of the vortex core structure should be possible because an 
advantage here is that the length scale of the core radius 
is much larger than those in superfluid "^He or supercon- 
ductors. 
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V. CONCLUSION AND DISCUSSIONS 

We have investigated the stability conditions of a vor- 
tex in BEC systems by performing microscopic self- 
consistent calculations based on the Popov approxima- 
tion. It is demonstrated that the BEC systems confined 
in a harmonic potential at finite temperatures sustain the 
stable vortex. This is contrasted with the cases in the ab- 
sence of the pinning potential where the same calculation 
yields unstable vortex as shown previously [|l^,^. The 
vortex instability comes from the negative eigenvalue of 
the LCLS (lowest core localized state). Therefore phys- 
ically it is possible to lift this energy by increasing the 
effective potential felt by this quasi-particle. This can be 
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FIG. 1. (a) The condensate density \4>{r)f at various T. 
The total number of the condensate |<?i(r)|'^ decreases as T in- 
creases. Inset shows the profile of \(j){r)\'^ at small r. (b) The 
non-condensate density p(r). It grows as T increases. The 
characteristic length scale near the center corresponds to the 
core radius of |(/>p in (a), (c) The total density 10(^)1^ + p{r). 
Each figure shows the densities at T = 200, 500, 700, and 
900nK. 
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FIG. 2. The dispersion relations of the eigenvalues Sq at 
T = 200nK. (a) Along qe. The black points are states with 
= 0. The lowest black point at each qe corresponds to 
qr ~ 1, and the second lowest to qr = 2, and so on. (b) The 
eigenvalues along q^. Only the qe — —1 states are plotted. 
The solid curve is explained in the text below the Eq. (po[). 
The dashed curve denotes ^gf- 
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FIG. 3. The wave functions u{r){qe = —l,qz = 0) at 
T — 200nK. The radial quantum number qr varies from 1 
to 5. Each u{r) has node(s) whose number is equal to qr- 
These states are localized at the core. The vertical axis is 
arbitrary unit. 
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FIG. 4. The effective potential Vcs{r) at temperatures 
T = 200, 500, 700, and 900nK. We subtract cfiemical poten- 
tiai n from Kff(r) of eacli T. 



c 

CO 











-5 






-3 




(-1) 


qe=i 


-1 







200 



400 



800 



1000 



600 
T[nK] 

FIG. 5. The eigenvalues e of the selected lowest modes 
as a function of T. The number of each line denotes qe- All 
of the lines are eigenvalues with (g^ = l,gz = 0), except the 
line labeled "(—1)" with q — (2,-1,0). Solid lines mean 
gs < 0, dotted lines ge > 0, and the dashed line gg = 0. 
While the dipole mode with gg = 1 is exceptionally almost 
T- independent and fixed to e = 1, coinciding with the trapped 
frequency, other modes with qe > {qe < 0) tend to increase 
(decrease) with increasing T. The eigenvalues with large qe 
which are not shown here are generally increase. This is be- 
cause these wavefunctions generally lie at large radius, and 
the effective potential VcB{r) there increases as increasing T. 
See the large radius side of Fig. H. 
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ico 

FIG. 7. The local density of states N{r, E) with the lower 
energy side (0 < £ < 7 in trap unit) at T = 2QQnK. The 
states distinctively localized near the core correspond to the 
angular momentum qe = —1. The sharp peak at the left (low 
energy) side is the LCLS. Some too high peaks around the 
core are cut in a certain value of N(r,E). 
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FIG. 8. The non-condensate contribution j^^\r) of the 
particle current density at T = 200nK. Each component with 
qe ~ 0, ±1,— 2 is shown separately. Inset shows the total 
current density je (r) 



(r)+jf' (r). 
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FIG. 6. T dependencies of ^'s in various definitions. Both 
^ (solid line) and ^peak (lower dashed line) are defined in the 
text. The upper dashed line means 1.5 x ^peak- The dotted 
line labeled is the peak radius of je{r). 
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FIG. 9. Various properties of a vortex system stabilized 
by the pinning potential in Eq. (H). Vo = 5 (trap unit), 
ro = 1/im, and T = OK. (a) The total particle number den- 
sity, (b) The dispersion relations of the eigenvalues e along 
qe- Black dots are eigenvalues with qz ~ 0. (c) The wave 
functions u{r){qr — 1 to 5, = —l,qz ~ 0). The vertical 
axis is arbitrary unit. 
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FIG. 10. Various properties of a vortex system stabilized 
by the pinning potential. Vo ~ 50 (trap unit), ro = l/im, 
and T = OK. (a) The total particle number density, (b) 
The dispersion relations of the eigenvalues e along qg. Black 
dots are eigenvalues of q^ — 0. (c) The wave functions 
u{r){qr = 1 to 5, = —l,qz = 0). The vertical axis is ar- 
bitrary unit. 
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